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Abstract 

The corrections of the gluon fusion to the BFKL equation in a unified 
partonic framework are studied. This modified BFKL equation predicts 
a stronger shadowing, which suppresses the gluon density and even leads 
to the gluon disappearance below the saturation region. We suggest that 
this unexpected effect is caused by a possible chaotic solution of the new 
equation. 
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1 Introduction 



The gluon density is an important knowledge in the researches of the high 
energy physics. The gluon distribution in the nucleon is mainly extracted from 
the quark distributions using the DGLAP equation [1,2] since the gluons are 
uncharged [3]. However, the above mentioned gluon distribution is unreliable 
at the small x region. The BFKL equation [4] should replace the DGLAP 
equation at small x, since the former sums up the leading powers of \og(l/x) 
(the LL(l/x) approximation). However, the rapid growth of the gluon density 
in the DGLAP and BFKL equations arises the gluon fusion, which leads to 
the nonlinear corrections in these equations. 

There are several nonlinear evolution equations considered the corrections 
of the gluon fusion to the DGLAP and BFKL evolutions. One of the most 
widely studied such models is the Balitsky-Kovchegov (BK) equation [5]. The 
BK equation is written for the scattering amplitude. It can be phenomeno- 
logical explained by the dipole model [6], in which the nonlinear terms are 
formed by the dipole splitting and the screening effects origin from the double 
scattering of the probe on the final states. The BK equation is equivalent to 
the leading part of the decoupled JIMWLK equation [7], which is an infinite 
hierarchy of the coupled evolution equations for the Wilson lines. A remark- 
able feature which emerges from the solution of the BK equation is that the 
scattering amplitude gradually approaches to a limit form. This behavior is 
called the saturation, where the gluon fusion balances with the gluon splitting. 
The researches to find the saturation effects has been very active in the past 
few years. 

On the other hand, the nonlinear dynamics may have a characteristic 
solution-chaos, which have been observed in many natural phenomena [8]. 
Therefore, it is interesting to ask: Does nonlinear QCD evolution equation of 
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gluon distributions also have chaotic solution? And how does it affect on gluon 
distributions? 

In this paper we present a nonlinear modified BFKL (MD-BFKL) equation 
in a unified partonic framework, and this equation describes the corrections of 
the gluon recombination to the BFKL equation but differ from the BK equa- 
tion. We find that the unintegrated gluon distribution function F(x, k 2 ) in the 
MD-BFKL equation begins its smooth evolution under suppression of gluon 
recombination likes the solution of the BK equation. But when x comes to 
a critical x c , F(x, k 2 ) will oscillate aperiodically near the evolution endpoint 
k 2 . This is the characteristic feature of Chaos: random and sensitivity to ini- 
tial condition. Furthermore, with the enhancement of oscillation, distribution 
F(x, k 2 ) will disappear suddenly at x c . 

Our idea is straightforward. We begin with an elementary amplitude Fig. 
la of deep inelastic scattering (DIS) in the parton scattering picture, where the 
dashed line implies a virtual current probing gluon. This amplitude, together 
with its conjugate amplitude, constructs the DGLAP equation for gluon. The 
correlations among the initial partons are neglected in the derivation of the 
DGLAP equation. Obviously, this assumption is invalid in the higher density 
region of partons, where the parton wave functions begin to spatially over- 
lap. Therefore, the corrections of the correlations among initial gluons to the 
elementary amplitude at small x should be considered. To this end, we add 
a possible initial gluon to Fig. la step by step. The resulting three sets of 
amplitudes are listed in Fig. lb-Id. It is interesting that these amplitudes 
produce the BFKL equation [4], the modified DGLAP equation [9] and a new 
BFKL (MD-BFKL) equation, respectively. 

As we well know physical picture of the same evolution process has different 
appearances in different coordinate systems. We shall expand our derivations 
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of the evolution equations in the Bjorken frame, where the traditional parton 
distributions inside a fast moving target are exactly defined in the factorization 
scheme. However, the BFKL equation was originally produced by using the 
ladder diagrams of reggeized gluons [4]. Although the BFKL equation can 
be reproduced in the dipole model in the target rest frame [6], we try to put 
the BFKL equation, the DGLAP equation and as well as their modifications 
together into a same partonic framework as shown in Fig. 1. Therefore, at 
first we reconsider the derivation of the BFKL equation in Sec. 2, where the 
time ordered perturbative theory (TOPT) [10] is used. 

The evolution kernel of the DGLAP equation has singularities (i.e., in- 
frared (IR) divergences), which relate to the emission or absorption of quanta 
with zero momentum. However, the singular terms provide the leading con- 
tributions to the DIS processes. Since a correct theory is IR safe, the IR 
divergences are cancelled by combining real- and virtual-soft gluon emissions, 
and the leading contributions are retained. Altarelli and Parisi first apply a 
1/(1 — z) + prescription to the singularities so that the integrals in the evo- 
lution equation are finite. In the same time, they add to each kernel a term 
proportional to 5(1 — z) to keep the momentum conservation [1]. The above 
mentioned treatment by hand was advanced by Collins and Qiu in a covariant 
field theory [1 1] . More simple calculations of the virtual diagrams are proposed 
via the TOPT in [9,12]. We find that the above mentioned real and virtual 
diagrams have a similar time ordered structure but with different cut lines. In 
general, we use a probe to observe the parton distributions inside the target 
in DIS, we cannot control the probing positions. In principle, we should sum 
over all cut diagrams belonging to the same time-ordered un-cut graph, and 
these graphs have a similar singular structure but may come up with oppo- 
site signs. On the other hand, the DIS structure functions are the imaginary 
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parts of the amplitudes for the forward 'Compton' scattering of the target 
with a probe. Using the time-ordered perturbative expansion of the statement 
of unitarity of the S-matrix, one can prove that the structure functions are 
associated with the sum of cut diagrams. These different cut graphs represent 
various possible sub-partonic processes due to the unitarity of the perturbative 
S-matrix. Therefore, the sum of cut graphs is necessary not only for infrared 
safety, but also for collecting the leading contributions and restoring unitar- 
ity. For this sake, the TOPT-cutting rules are proposed to present the simple 
connections among the relating cut-diagrams including the real- and virtual- 
diagrams in our previous works. We shall detail this method in Appendix A. 
The BFKL-kernel also has singularities on the transverse momentum /c-plane. 
Thus, we can pick up the contributions from the virtual diagrams using the 
TOPT-cutting rules without the complicated calculations in the re-derivation 
of the BFKL equation. 

The new derivation of the BFKL equation allows us to add the corrections 
of the gluon fusion on it successively according to the physical pictures in 
Fig. 1. Such derivation of the MD-BFKL equation is given in Sec. 3. In 
the nonlinear part of the MD-BFKL equation there is similar IR divergences 
as in the BFKL-kernel. Naturally, the same regularization scheme as in the 
BFKL equation is necessary, i.e., we use the TOPT-cutting rules to pick up 
the contributions from the virtual processes in the linear and nonlinear parts 
of the MD-BFKL equation. This leads to the MD-BFKL equation, which is 
different from the BK equation. 

The discussions of the MD-BFKL equation are shown in Sec. 4. We find 
that the MD-BFKL equation has the close relationship with the DGLAP, 
BFKL and MD-DGLAP equations. We shall present it in Appendix B. The 
differences between the MD-BFKL equation and the BK equation are empha- 
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sized. We point out that the evolution kernels in the linear and nonlinear 
parts of the MD-BFKL equation are regularized by using the same method: 
the summations of the real and virtual processes. While these singularities in 
the linear and nonlinear parts of the BK equation are regularized by using two 
ways: the former is via the contributions of the virtual diagrams, the latter 
absorbs these singularities into a definition of the double amplitude. This leads 
to the different predictions of the two evolution equations. The numerical so- 
lutions of the MD-BFKL equation show a unexpected solution: the shadowing 
effect suddenly increases near a critical small value x c . This stronger shadow- 
ing breaks the balance between the gluon fusion and splitting. In consequence, 
the gluons disappear at x < x c below the saturation region. It is interesting 
that some properties of the solutions impels us to suggest that this unexpected 
effect is caused by a possible chaotic solution of the MD-BFKL equation. We 
shall detail them in Sec. 4. Obviously, the present form of the MD-BFKL 
equation is derived in a partonic framework, which contains some unclear con- 
ditions. As an available QCD evolution equation, the MD-BFKL equation is 
needed to further justify. We shall discuss them in the last section. 



6 



2 The BFKL equation 

We consider following partonic picture for the DIS process. At the lowest 
order, the elementary amplitude in Fig. la together with its conjugate ampli- 
tude constructs the DGLAP equation for gluon. However, this picture should 
be modified at small x due to the correlations among initial gluons. For exam- 
ple, the lowest corrections to the DGLAP-amplitudes are given in Fig. lb, or 
their details in Fig. 2, where the dashed lines in Fig. 2 are the time-ordering 
lines in the TOPT. These processes imply a scatted gluon before its radia- 
tion is omitted from two correlating gluons. We call such a correlating gluon 
cluster as the cold spot, which phenomeno logically describes the correlation 
among initial partons. The elementary correlations among two initial partons 
in a cold spot are shown in Fig. 3, where the dark circle implies all possible 
Q CD- interactions. 

The evolution kernel in a QCD evolution equation is a part of more compli- 
cated scattering diagram. In general, the correlations among the initial partons 
make these partons off mass-shell. One of the most important hypotheses in 
the derivation of the DGLAP equation is that all correlations among initial 
partons are negligible during the interaction. Therefore, the interaction of a 
virtual probe with the nucleon can be factorized as the nonperturbative par- 
ton distribution and hard probe-parton matrix in the collinear factorization 
scheme. At the higher density range of partons, these correlations among the 
initial partons can no longer be neglected. In this case, the transverse mo- 
menta of the initial partons are non-zero and these partons are generally off 
mass-shell, therefore, the fey-factorization scheme is necessary. 

In this work we use the semi-classical Weizsacker- Williams (W — W) ap- 
proximation [13] to realize the /^-factorization scheme. The reason is that the 
W — W approximation allows us to extract the evolution kernels and to keep 
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all initial and final partons of the evolution kernels on their mass-shell. For 
illustrating this idea, we consider a sub-process as shown in Fig. 4, where the 
momenta of three massless particles are parameterized as 



Pi = (p,Q,p), 



k 2 

P2 = (x 2 p+ 7 ^,t,x 2 p), 
2x 2 p 



and 



Ps = (p°3,-k,x 3 p), (2.1) 

where p — > oo. The particle 3 is off its mass-shell. We decompose this propa- 
gator into a forward and a backward propagator according to the TOPT. The 
corresponding momenta are on mass-shell 



k 2 

Pf = (x 3 p + -k, x 3 p), 

2x 3 p 

and 



k 2 

Pb = (-xap - tt^, -k, x 3 p), (2.2) 
2x 3 p 

respectively. Note that Eqs. (2.1) and (2.2) require 

xp^>\k\, (2.3) 

which is consistent with the conditions of high energy- and small angle-scattering 
in the W — W approximation. One can find that the forward propagator con- 
tributes the power-factor 1/fc 2 , since 



111 . . 

(2.4) 



2E 3 E 2 + Ep F - E x k 
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while the backward propagator gives the contributions ~ 1/p 2 - If the domi- 
nant contributions in a process are from the terms with the power of transverse 
momentum, one can find that the contributions from the backward propaga- 
tor are suppressed at p — > oo, and the dominant contribution is from the 
forward component. In this case, we can break the forward propagator since 
it propagates on mass-shell. 

According to the scale-invariant parton picture of the renormalization group 
[14], the observed wave function ^(x 2 ,k) is evolved from the initial wave func- 
tions ty(xi,p a ) and ty(xi,p b ) via the QCD interactions, i.e., 





*(x 2 , k) = ^(x 1 ,p a )A BFKL1 + *(xi,p 6 )A BFKL2 , 



(2.5) 



where two perturbative amplitudes corresponding to Fig. 2 are 




(2.6) 



and 




(2.7) 



The momenta of the partons are parameterized as 




(2.8) 



k = (x 2 P + 



k 2 



,k, x 2 P), 



(2.9) 



2x 2 P 



l a = ((x 1 -x 2 )P + 



2( Xl -x 2 )P 




,l_ a ,(x 1 -x 2 )P), 



(2.10) 




(2.11) 
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and 



h = ((xi - x 2 )P + - Lb k, ( Xl - x 2 )P). (2.12) 

Z[Xi — X 2 )r 

The matrixes of the local QCD interactions are 



Mi = igf abc [g a p{p a + k) 1 + g Pl {-k + l a ) a + g 7a (-l a - P^pYaipa^pify^iQ, 

(2.13) 

M 2 = igf abc [g a p{pb + k) 1 + g Pl {-k + l b ) a + g ia {-l b - Pb)^a(Pb)^(k)e 7 (l b ), 

(2.14) 

where the polarization vectors are 



e(p a ) = (0,e_,- -' { ^ + p L) ), (2.15) 



C (fc) = (0, £ ,-^|), (2.16) 

and 

e(O = (0,e,- ? ), (2.17) 

where e is the transverse polarization of the gluon in e M = (eo, e, e^j = (0, e, 0), 
since the sum includes only physical transverse gluon states in the TOPT form 
[1,10]. 

Taking the LL{l/x) approximation, i.e., assuming that x 2 C xi, one can 
get two similar amplitudes 



A BFKLl = tgf abc 2M^, (2.18) 
V x 2 k 
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and 



A B fkl2 = igf abc 2M^-J=. (2.19) 
\l x 2 k 

However, these two amplitudes really occupy different transverse configura- 
tions. This is a reason why the dipole model of the BFKL equation is derived 
by using the transverse coordinator-space. However, we shall show that the 
momentum representation still can be used to distinguish the differences be- 
tween Eqs. (2.18) and (2.19). 

The two parton correlation function is generally defined as 



I 1 2 = f( x ,Pa,P b ) 

= f^Ea±Ek iSa -p^ = f( X ,k c ,k ab ), (2.20) 

where fc c and are conjugate to the impact parameter and transverse scale 
of a cold spot. Equation (2.20) implies the probability finding a gluon, which 
carries the longitudinal momentum fraction a; of a nucleon and locals inside a 
cold spot characterizing by k^ and /c afe . 

In this work we derive the evolution equations in the impact parameter- 
independent case. This approximation implies that the evolution dynamics of 
the partons are dominated by the internal structure of the cold spot. Thus, 
the evolution kernel is irrelevant to /c c and we shall use 

f(x,k ab )= (2-21) 
which has the following TOPT-structure 



E. 



ah 



2E f 



M. 



1 ^ ' 



P-*k ah X 



_Ep — E, 



ab 



Ey 



IE, 



ab- 



n 



d 3 k 



x 



(2nf2E x - 



[2.22) 
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Notice that all transverse momenta in Eqs. (2.5)-(2.19) are indicated rela- 
tive to the mass-center of the nucleon target. However according to Eq. (2.21), 
the evolution variable is the relative momentum A^ b , therefore, it is suitable 
to rewrite all momenta to relative to p b in Eq. (2.6) and to p a in Eq. (2.7), 
respectively. Thus, we replace the transverse momenta as follows: 

k ~^ k~ = kobi 



and 



in Eq. (2.6) since 



La ~^ kab kob — P n k — kaQ, 



kab — kaO + kddb, 



(2.23) 



(2.24) 



and 



and 



k ~^ k~ p n = koai 



Lb ~~ ^ kba koa — P b k~ khOi 

in Eq. (2.7). In consequence, we have 



(2.25) 



and 



(2.26) 



A B FKL(kaO:kob:Xl,X 2 ) = igf abC 2J — 

V x 2 



k 2 



(2.27) 
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where we identify two e in Eq. (2.27) since the measurements on (x 2 , Mo) an d 
(^2,&ofe) are really the same event. 

Equation (2.5) with Eqs. (2.26) and (2.27) provide such a picture: a parent 
cold spot with the longitudinal momentum fraction X\ and transverse momen- 
tum k^ b radiates a gluon, which has the longitudinal momentum fraction X2 
and the transverse momentum k^Q (or /c ofl ). It is interesting that this is a 
picture like the dipole model but in the full momentum space. In fact, using 
the Fourier transformation, one can obtain the corresponding amplitude in the 
dipole model [6] 



A B FKL(x a0 ,x 0b ,x 1 ,x 2 ) = J ^Tr) ^ AbFKL ^- a0 ' - 0? " Xl ' %2 ) 6 ~°° ^°° + iM ° b ~°" 

= ^-2 S|f + |^].e... (2.28) 

where x is the conjugate coordinator corresponding to the relative transverse 
momentum k. 

We take the square of the total amplitude, one can get 



«6 



2Er 



Mp^ Kh x\ 



1 



---p — E a b — E x 



1 



2£. 



ab 



2 



n 



(2tt) 3 2£ x 



/I A* a " 

a bfklA bfkl - 



pol 



ab 



1 d 2 k' 
x— — iM^fc^f (2vr) 4 5 4 (g profee + A; - A/) 



X\ 

= f(xi,kab) ® —K-BFKL (tab,kiO,®s) ® d<j(q* robe k(x 2 , k^) -> tf(x 2 ,U)) 
X2 

= A[^(x 2 ,k a0 )^*(x2,k a0 ) + *(x 2 ,^ )**(^2, £ &)+ 
^(^2,io?,)^*(^2,iao) + *(z 2 ,i^&)**(~-2,iy)&)] ® da(q* robe k(x 2 , hao) -»• fc'(x2,fc')) 
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= Af(x 2 ,k a0 ) <g> da(q* robe k(x 2 ,t a0 ) -> fc'(ar 2 , £'))> ( 2 - 29 ) 

where on the last step, because the probe only picks up the contributions from 
^(%2, kao)^*( x 2, kajo), we regard Af(x 2 ,kao) as the increment of the distribu- 
tion f(xi,kab) when it evolves from (xi^kab) to (x 2 ,k a0 ). Therefore we have 



Af(x 2 ,k a0 ) 

= I d jT- I ——^BFKL(k ab ,k a0 ,as)f(x 1 ,k ab ), (2.30) 

J kab Jx 2 x l x 2 

or 



AF(x 2 ,k a0 ) = Ax 2 f(x 2 ,k a0 ) 
-jT 1 / JC B FKL{kcb,k 1 o,ois)F(xi,k^ ) ). (2.31) 

tab Jx 2 X l 

Using definition 

F{x 2 ^ Q ) = F(x 1: k^ b ) + AF(x 2i y, (2.32) 

we write [9] 

dx 

= J d^ lb ICBFKL(kab,kaO,as)F(x,k ab ), (2.33) 

According to Eq. (2.26), the evolution kernel reads as 

%r n i . x \dx\ 2 V* /I a* dx\ d 2 kab 

>^BFKL{Kab,K<iO, a s) « Kab ~ 2^ A BFKlA bfkl ~ 

^2 x l po i * x l \*Tt) 

Finally Eq. (2.33) becomes 
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x dF(x,k a0 ) 
dx 

= ^ f^j&rFix^). (2.35) 

This is the real part of the BFKL equation. 

Using the TOPT-cutting rules [9,12], one can prove that the virtual dia- 
grams in Fig 5 contribute the similar evolution kernel as the real kernel but 
differ by a factor —1/2 x (1/2 + 1/2) (see Appendix A). Here we call the cut 
diagram with a naive partonic amplitude without any QCD-corrections as the 
virtual diagram. The negative sign arises from the changes of time order in 
the energy denominators. The factor (1/2 + 1/2) is due to the fact that the 
probe "sees" only the square root of the parton distribution accepting the con- 
tributions of the partonic processes in a virtual diagram, and the other factor 
1/2 origins from the symmetry of the pure gluon process [9]. Therefore, the 
evolution equation corresponding to Fig. 5 is 



dF(x,k ab ) 
-x- 



i ^0 ,2 , v Hx,ko>). (2.36) 



dx 

Since we calculate the contributions to AF(x, ka ), we should make the re- 
placement b <-> in Eq. (2.35). Combining the real and virtual parts of the 
evolution equation, we have 



a s N c 



x dF(x,k a0 ) 
dx 

k 2 ~ k 2 
-ab 77/ i \ ^a0 



F(x,k ab ) — ^° 2 F(x,k a0 ) 



k 2 k 2 v '- aby u 2 h 2 

t±ab±«b 



(2.37) 



According to Eq. (2.22), the distribution f(x,k) in the TOPT-form con- 
tains a singular factor 1/k 4 , which arises from the off energy-shell effect in the 



15 



square of the energy denominator. In order to ensure the safety of taking the 
W — W approximation [13], we remove this factor to the evolution kernel and 
use the following new definition of the unintegrated gluon distribution 



k 4 ~ 

F(x,k) = ^F{x,k), (2.38) 
k 

where k is a unity vector on the transverse momentum space. Thus, Eq. (2.37) 
becomes 

x dF(x,k a0 ) 
dx 

= ¥t[ *k*Jkr [2F(x, U - F(x, fa)] , (2.39) 

which consists with a usual form of the BFKL equation. 

The correlations among the initial gluons can be neglected in the dilute 
parton system. In this case the contributions of the interference diagrams 
Figs. 2c and 2d disappear. Thus, the kernel (2.34) reduces to the splitting 
functions in the DGLAP equation at the small x limit, 

xi dxi 2 a s N c dxidk 2 

I^BFKL{Kab, KaO, a s) « Hab ~> TY 

X 2 Xi 7T X 2 k 

= fcDGLAp-F^ -■ (2.40) 

k xi 

The two initial gluons have the same transverse momentum and we always 
can take it to zero and use the collinear factorization to separate the gluon 
distribution. The corresponding DGLAP equation reads [12] 

2 <9#(x B ,Q 2 ) f 1 dx, fx B \ 2 

Q dQ 2 = I — ' Cdclap W' a 'J g(Xl ' Q > 

71 Jx B X 1 X B 
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where the scaling restriction 5(x 2 — xb) is included and 



rQ 2 dk 2 rQ 2 dk 2 

G(x,Q 2 ) =xg(x,Q 2 ) = I ^xf(x,k 2 )^l ^F(x,k 2 ) (2.42) 

is defined. 

The above derivations reveal a relation between two evolution equations 
from a new view point: the evolution kernel of the DGLAP equation is the un- 
interferant part of the BFKL kernel. Or reversely, the BFKL equation can be 
regarded as the leading corrections from the correlations of the initial partons 
to the DGLAP equation. 
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3 The modified BFKL equation 

We consider the evolution kernel based on Fig. Id, which constructs the 
MD-BFKL equation. Notice that two pairs of initial gluons, which are hidden 
in the correlation function, for example in Fig. 6a, should be indicated as Fig. 
6b. A set of cut diagrams based on Fig. Id are listed in Fig. 7, where the 
probe vertices have been separated out using the W — W approximation [13]. 
Similar to the derivation of Eq. (2.34), we write the evolution kernel of the 
MD-BFKL equation as 

K-MD-BFKL — ' „ 2 A-MD— BFKL A*MD— BFKL- (3-1) 

Xl pol 

The amplitudes 



AmD-BFKL — AmD-BFKLI + AmD-BFKL2, (3.2) 

where 



Amd-bfkli - — Zr^TTTT p 7 p w ^-M MD _ BFKL1 , (3.3) 



and 



/ 2,E 1 1 
A M d-bfkl2 = J k 2E,E, + E, -E -E Mmd ~ BFKL2 - ^ 

The momenta of the partons, for example, are parameterized as 



(I -m) 2 

Pa = ( Xl P + 2a; ~ ,L ~ Uk, xiP) , (3.5) 



p b = {xiP + ( ~ 2 ^p )2 ,fe + m, xiP), (3.6) 
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k 2 

k = (x 2 P+^-p,k,x 2 P), (3.7) 



l a = {{2 Xl - x 2 )P + 2{2xi L ^ X2)p ,L, (2xi - x 2 )P). (3.8) 



Pc = (xiP + ( ~ + ^ )2 ,fc + m', xiP), (3.9) 



p d = (x 1 P+ K - a ^-m'^P), (3.10) 

IX\P 



l d = {{2 Xl - x 2 )P + 2{2x ^_ X2 ) P ^ ^ ~ X2 ) P )' ( 3 - U ) 
For example, in the t-channel 



(]c -\- m) 2 k 2 

m = p b -k = ((xi - x 2 )P + K ~ 2x ~ J - 2 T x z ^,m, (xi - x 2 )P), (3.12) 

and 



m' = p c - k = ((xx - x 2 )P + ( - 2 + 2 g° 2 - m', (x! - x 2 )P). (3.13) 

The matrixes in Eqs. (3.3) and (3.4) are 



-iff* 

M MD . BFKL1 =^^/^C'^— ^^^/^C^e^pJe^p^e^/Je^A;), (3.14) 
and 



-id 1 * 1 

M MD _ BFKL2 = igf ABC C a ^^^igf^C^e a { Pd )e p {p c )e%l d )el{k), (3.15) 
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where d?^ = rVrf 1 + rfn 1 — g lv , C a ^ 1 C pa ' n are the triple gluon vertices and 
polarization vectors are 



e(Pa) = (0,e,- £,( ^ p m) ), (3.16) 



e( Pb ) = (0,e,- € -- {k+ p m) ), (3.17) 



e(fc) = (0,e,-^|), (3.18) 



and 



Thus, at small x we have 

A M D-BFKL(k, Xl, X 2 ) 



= 9 2 f ABC f DCE J^- 2 



^ € ' k € ' fc 6 ' Jv € ' Jv 

6_r 2 7T~ + 6_ T2 TT" 



fc fc k 



(3.20) 



where one of two factors in every term is from the approximation 

e(k)m/m 2 ~ e • A;/fc 2 , 

and 

e(k)m'/m' 2 -e-k/k 2 . (3.21) 

In the TOPT calculations, the contributions of the propagator with the 
4-momentum m (or m') can be separated from that of the propagator with 
the 4-momentum k, and those propagators belong to the different cold spots, 
respectively. We use the relative transverse momenta to replace the relating 
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momenta in Eqs. (3.5)-(3.13) and recalculate Eqs. (3.3) and (3.4). The result 
is 



9 2 f ABC f DCE J^- 2 



A-MD—BFKL (it) Xi, X 2 ) 

e ■ k/ 



- ' —(Pb,Pc)or(p a ,p d ) 



» e-k (pbiPc) 



+6 



e • kr 



£ • Z-(p b ,Pc)or(p a ,p d )( m ') e • k 



-(Pb,Pc)or(p a ,p d )( m ) l±(p b ,Pc) 
£ ' k-(p b , Pc ) 
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(3.22) 



where the foot-indexes of the relative transverse momenta indicate the corre- 
sponding cold spots and k(m), k(m') imply that the momenta origin from to, 
m', respectively. Using the definitions 



— bc ~ Eb ~~ Re' — 60 ~~ Rh~ — > — oc — k ~ P r -> 



(3.23) 



we have 



kbc — kbO + hoc- 



(3.24) 



We read two momenta k( PbtPc ) in Eq. (3.22) as k 0c and k b0 , respectively. On 
the other hand, due to momentum conservation, we have 



-(Pb,Pc) 



m 



P h -k = Pl l ~Pr--k + Pr- = kbc-koc = kbO> 



a> tic 



and 



hpb, P c) M =k-E. = k-p h -p J . + p h = k b-kcb = hoc- 



Thus, we obtain 



(3.25) 
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V ^ X2 ^feO-^Oc 

Note that two factors kl and kl c in the denominator of Eq. (3.26) are corre- 
lated through Eq. (3.24) and they have double poles as in the BFKL-kernel 
(2.34). 

The result Eq. (3.26) seems irrelevant to p a and p d . However, there are 
two possible contributions of the cold spot (p a ,Pd) to the evolution kernel: 

(1) The momenta p a and p d don't flow into the amplitude (3.26). Therefore, 
the cold spot (p a ,Pd) in Fig. 6b is independent of the evolution dynamics and 
its distribution should be integrated as a unobservable quality. Thus, the 
resulting kernel reduces to the linear BFKL kernel. 

(2) The momenta p a and p d flow into the amplitude (3.26) through m and 
rnf. The momenta k( PbtPc )(m) and k( PbtPc )(m') in Eq. (3.22) are alternatively 
replaced by k (PajPd) (m) = p a - k = fc a0 and k {pa>Pd) (m') = k - p d = k od , 
respectively. The corresponding amplitudes become 



A' 

^MD-BFKL 



a 2 fABC fDCE j x l <L' ^aO^ " kp c (Q Q7\ 



and 



A" 

^MD-BFKL 



= Qg 2 f ABC f DCE jK ^- k f-- 2 kb \ (3.28) 

where one can introduce 

^ab = P a -P d = P a -k-p d + k= k aQ + kjoa. (3.29) 

In general, the momenta and /c M in Eqs. (3.27) and (3.28) are undetermined 
since l a and l d in Fig. 6b are unobserved, they should be integrated out as 
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two independent variables. Thus, the resulting evolution kernel reduce to the 
DGLAP-like kernel. 

Obviously, the above mentioned two situations should be excluded in our 
resummation for getting the leading corrections, unless we have following re- 
striction conditions 



kod — hoc-, (3.30) 



and they imply that 



kad — kba (3.31) 

due to Eqs. (3.24) and (3.29). For understanding Eq. (3.31), we image that 
before the probe interacts with the target, two overlapping cold spots have 
recombined into a common cold spot (p b ,P ), i.e., Fig. 6(b) — >■ Fig. 7(c). This is 
an inverse processes of the dipole splitting in the BK equation [4,5]. Therefore 
the probe always measures the recombination processes of four initial gluons, 
which origin from a same cold spot and share a same relative momentum. 

Summing all channels, we derived the evolution kernel corresponding to 
Fig. 7 and the result reads 



y~MD—BFKL dX\ j2j z 



A M d-bfklA* md _ bfkl 
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X 2 Xi 

1 t], \n 



1 dx\ 2 



18a 2 s N 2 C 1 t£ c dxi 

vr iv c 2 -ig;gs^ 



167T 3 X\ 

; "'I'hr- (3.32) 
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In the case decreasing gluon density, the contributions of the interference 
terms (Figs. 7c and 7d) disappear and Fig. Id return to Fig. lc. Thus, Eq. 
(3.32) reduces to the real part of the MD-DGLAP- kernel [9] 



K MU - a "^——d z k 



xi dxi 2 18a^ N 2 dx\d 2 k 

x 2 x\ 7i N; — 1 x 2 k 



— y~MD—DGLAP dX\ dk ^ ^ 

X\ k_ 



Similar to Eq. (2.33) we derive the MD-DGLAP equation using 



G(x 2 ,Qj) = G(x u Ql) + AG(x 2 ,Q 2 2 ) 

= G(x 1 ,Q 1 )+ — j- / ^MD-DGLAP ( j 0* S )G (> (x±, Qi), 

J Qlmin Ql -W 2 Xi Xi \Xi J 

(3.34) 

where a power suppressed factor l/Q\ has been extracted from the evolution 
kernel. 

The 4-gluon correlation function G^> is a generalization of the gluon dis- 
tribution beyond the leading twist. It is usually modeled as the square of the 
gluon distribution [9]. For example, 



G^ 2 \x,Q 2 ) = ^G 2 (x,Q 2 ), (3.35) 
where Rn is the effective radius of nucleon. 



The complete MD-DGLAP equation includes the contributions of the two- 
parton-to-two-parton (2 — > 2) amplitudes, the interference amplitudes between 
the one-parton-to-two-parton (1 — > 2) and the three-parton-to-two-parton 
(3 — > 2) amplitudes. Where we meet very complicated calculations about the 
interference- and virtual amplitudes. However, the TOPT-cutting rules show 
that the above mentioned amplitudes correspond to the same recombination 
kernel while they work on the different kinematic regions [9]. In particular, 
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the contributions of the virtual diagrams are cancelled with each other in the 
MD-DGLAP equation. 

Other key problem is that we meet various multi-gluon correlation func- 
tions, in which the cut line cuts off the nonperturbative matrix with different 
ways. Fortunately, Jaffe has shown that these correlation functions on the 
light cone has the same form in the DIS processes [15]. The Jaffe-cutting rule 
was broadly used in the study of the high twist processes. The TOPT provides 
a straightforward explanation about the Jaffe-cutting rule: since all backward 
propagators are absorbed into the nonperturbative correlation functions in the 
collinear factorization scheme, the partons correlating two initial gluons inside 
the nonperturbative matrix are on mass-shell. Therefore, the correlation func- 
tions with cuts at different places are the same [16]. Thus, the Jaffe-cutting 
rule can be included in our TOPT-cutting rules (see Appendix A). Combining 
DGLAP dynamics at small x, the MD-DGLAP equation reads 

dG(x B ,Q 2 ) 



dlnQ 2 

7T Jx B Xi HR N Q 2 N 2 - 1 Jx B /2 Xi 

36a 2 N 2 rV*d Xl 2 2 
~ vR 2 N Q 2 N? - 1 L ~~x~[ (3 - 36) 

The second term on the right hand-side of Eq.(3.36) is positive antishadowing, 

while the third term is negative shadowing. 

Returning to the MD-BFKL equation. According to Eq. (3.31), the probe 

always measures the recombination processes of four initial gluons in a same 

cold spot. Thus, we model the 4-gluon correlation function as the square 

of the gluon distribution like in the leading twist case Eq. (3.35), i.e., 
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F (2 \x,k bc ) = J d 2 kR F (k bc ,k)F(x,k bc )F(x,k) = ^F 2 (x,k bc ). (3.37) 
Using the evolution kernel (3.32), we write 

F(x 2 ,k b0 ) = F(x 1} k bc ) + AF(x 2 ,k b0 ) 
~, , N 18a? f f i/2 d Xl 1 fc? — o . . 

= n*i,u + ^jjft-J^J^ -^wM F(X1 ^ 1 (3 - 38) 

Now let us discuss the contributions from the virtual diagrams. The TOPT- 
cutting rules show that the diagrams in Fig. 8 have a similar evolution kernel 
as that in Fig. 7 but with the different kinematical variables and differ from a 
simple numerical factor. 

The processes in Figs. 7 and 8 contributes the net positive antishadow- 
ing effect. The negative shadowing effect really origins from the interference 
processes, two of them are shown in Fig. 9. Here the contributions from the 
corresponding virtual processes also are necessary (see Fig. 10). The TOPT- 
cutting rules show that the processes in Figs. 9 and 10 also have a similar 
evolution kernel. 

Up to now we have separately established the relations of the evolution 
kernels between the real and virtual diagrams in the four to four partons 
amplitude and three to five partons amplitude, respectively. In the next step 
we show that a relationship between the above mentioned two kinds of virtual 
diagrams will link up all four evolution kernels. According to Eq. (3.20), the 
resulting amplitudes are irrelevant with the transverse momenta of the initial 
gluons at x 2 < ii. Thus, we use the relations showing in Fig. 11, which are 
derived in the collinear factorization schema [9] to reveal that the two kinds of 
virtual diagrams differ only from a minus sign, which is from a energy deficit 
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between two dashed lines in Fig. 10: because both the momenta A; M an d k 0c 
are indicated by k in the mass-center of the nucleon target, we have 



k 2 k 2 

1 >0, (3.39) 



2x m P 2x x P 

on the left hand-side of Fig. 10, where x m < xi, (x m and xi are the longitudinal 
momentum fractions in the momenta m and /, respectively); and 



k 2 k 2 



< 0, (3.40) 



2x m P 2 Xl P 
on the right hand-side of Fig. 11, where x m > x\. 

In consequence, we finally link up all evolution kernels and obtain following 
equation 

F(x 2 ,k m ) = F(xi,k bc ) + AF(x 2 ,k b0 ) 
fv , \ 18a s 2 N 2 r j2l rWdx! 1 k\ c ~ 2 . , . 

= F( Xl ,k ic ) + * C / d% c / — TTlfp-^ 0^1, fcftc 
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2 r 1 / 2 dx 1 1 k_l c ~ 2 

a !Hc , 2 ,2 ,2 r \ X lilhc) 



/x2 a?i k bc kj bQ k Jdc 
^r^^T^F^ho), (3.41) 

where we assume that the Jaffe-cutting rule still is held in the /^-factorization 
scheme (Fig. 12). The reasons are follows: (a) the propagator inside the cold 
spot is forward and on mass-shell at the W — W approximation; (b) the cor- 
relations to the cold spot from the other part of the unperturbative matrix 
are neglected in our model Eq. (2.21). Thus, the correlating function can 
be cut and we can use the same correlation function in the real, virtual, and 
interference processes. From Eq. (3.41) we have 
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dF{x,k ya ) 
dx 



18aj N 2 f 2 1 kj c ~ 2 fx 

9«s ?2/^ U \ f j2;„ 1 
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n 2 R 2 N N 2 -l v '- &U7 ^- 6c A:L^o 
Similar to Eq. (2.40) we note that 



F {2 \x,k hc )^ 
Thus we redefine 



E P — 2E bc — E 
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/V 4 



(3.43) 



F {2 \x,k) = J d 2 gR F {k,g)F{x,k)F{x,k!) 

= ^F^(x,k), (3.44) 
it 

where .Rp = -Rp/c jk . Submitting this equation with Eq. (2.37) to Eq. (3.42), 
the result is 



dx 

18«« N c [,2, 1 Ho F 2(X , 



V N 2 v2 /x , \ y^, 1 A; 2 
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Combining the linear BFKL equation, we finally obtain a complete MD- 
BFKL equation at small x 



dF{x,k m ) 
dx 

a*N r . f , , kin , s a s N r „, , , f , , 



^ 2 Y rf 2 fc fcc ^2-^2F(a;, - -^F{x,k J]0 ) J d 2 ^-^ , 



7T^ 2 JV? - 1 J Hi klk* * W - bC ) Ml W-l W J tin klk" 



(3.46) 



Comparing with the MD-DGLAP equation (3.36), the contributions of the 
virtual diagrams can't be cancelled in Eq. (3.46) and they are necessary for 
IR safety of the evolution equation. 
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4 Discussions 

We discuss the properties of the MD-BFKL equation. Like all the other 
models for the evolution equations, the derivations of the MD-BFKL equation 
also contain some assumptions. However, the rationality of the MD-BFKL 
equation can be immediately checked through its relations with other QCD 
evolution equations. For this sake, in Appendix B we summary a set of the 
evolution equations, which appear in this work. One can find that the MD- 
BFKL equation has the close relationship with the DGLAP, BFKL and MD- 
DGLAP equations. 

We don't calculate the MD-BFKL equation with the antishadowing correc- 
tions Eq. (3.46), which shall be studied elsewhere. For simplicity, we work near 
the saturation region, where the antishadowing terms in Eq. (3.46) disappear 
due to 

F(x,fc)~F(|,fc). (4.1) 
Thus, the MD-BFKL equation reduces to 

dF(x,k b0 ) 



-x- 



dx 



Oi<sN r f ,0, „„, , ,, cHf,N r , s f , , kj^ 



2 _ 2 j d 2 fc 6c -^^2F(x, fcfc) - -^F(x,k M ) J d 2 ]^ c -^r 



(4.2) 

An important difference between the MD-BFKL equation (4.2) and BK 
equation, which is usually written by using the scattering amplitude N(x, x) 
in the transverse coordinator space 



x dN{x m ,x) 
dx 
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= I d 2 x c X 2 m 2 [N(x bc , x) + N(x c0 , x) - N(x b0 , x) 

-N(x bc ,x)N(x c0 ,x)}, (4.3) 

is that they take the different regularization schemes: the singularities in the 
nonlinear real part are cancelled by the contributions from the corresponding 
virtual processes in Eq. (4.2), while such singularities are absorbed into the 
double amplitude N(x_ bc ,x)N(x c0 ,x) in Eq. (4.3). The following numerical 
calculations will show that the MD-BFKL equation has the chaotic solution, 
which suddenly arises a strong shadowing effect and cancels the gluons at small 
x. 

In our previous work [17], as an approximation, we used a known kernel 
K-md-dglap of the modified DGLAP equation to replace the kernel KLmd-bfkl 
in the MD-BFKL equation. In spite of this approximated approach, the result- 
ing equation nicely reproduces some characters of the BK evolution equation in 
the numerical results. For the convenience of following statement, we call this 
equation as the BK-like equation. Obviously, the above mentioned method is 
an approximation since in the MD-DGLAP approach there is no systematic 
resummation of small- x effects. In this article we shall derive the MD-BFKL 
equation, where the contributions of the gluon fusion at the LL{l/x) approx- 
imation are included. 

At first, we try to compare the MD-BFKL equation with the BK-like equa- 
tion [17] 



dF(x,k b0 ) 
dx 
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For the cylindrically symmetric solution we rewrite these two equations as 



dF(x,k 2 ) 
dx 

3a s k 2 [°° dti 2 \ F{x, y 2 ) - F(x, k 2 ) ^ F(x, k 2 ) 



81 a 2 r<*> dk' 2 k 2 F 2 {x, k' 2 ) - k' 2 F 2 {x, k 2 ) F 2 {x, k 2 ) . 



4 7cR 2 N Jg g 2 \ k! 2 W -k 2 \ ^4 ~ 4 ^/4 

and 



dF(x,k 2 ) 
dx 

,2 ,^ J7/2 f 77V „ l'2\ 171/ „ l2\ 771/ ,2\ 



_ aa s r r 00 dk I ^ (x, fc ) - b \x, k ) * {x,k ) I 81 a; 1 2 2 

(4.6) 

We use the Runge-Kutta method to compute the evolution equations. For 
example, we take the input distribution as 



F(x = I0~,k 2 ) = (3^k 2 /lGeV 2 exp 



\og 2 (k 2 /lGeV 2 )) 



(4.7) 



40 

where we take the normalized parameter f3 — 0.1. A similar input was used in 
[17,18]. In this work we use a fixed value of a s . Therefore, the cut-off in the 
integration can take a smaller value, for example, k£ = O.OlGeV 2 . 

The x-dependence of F(x,k 2 ) with fixed value of k 2 using the MD-BFKL 
equation (4.5) is illustrated by the solid curves in Fig. 13. For comparison, 
we draw the corresponding solutions of the BFKL equation via dashed curves. 
The results show that F(x,k 2 ) suddenly drops near a critical value of x ~ 
x c ~ 7.2 x 10~ 7 and this value is insensitive to the value of /c 2 . In particulary, 
the gluons disappear at x < x c . We try to explain the reason of the gluon 
disappearance as follows: 
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(i) We separately plot the contributions from the nonlinear shadowing 
terms of Eq. (4.5) in Fig. 14 and compare them with the results of Eq. 
(4.6) in Fig. 15. We find that the shadowing effect increases suddenly in the 
MD-BFKL equation near x C) while the shadowing effect grows gradually in Eq. 
(4.6). That is, the gluon disappearance is arisen by an anomalous shadowing 
effect in the MD-BFKL equation. 

(ii) To understand this strong shadowing better we compute the /^-dependence 
of F(x : k 2 )/k 2 with different values of x using Eq. (4.5) and illustrate the re- 
sults in Fig. 16. The dashed curves are the corresponding solutions of the 
BFKL equation. It is interest that the solution F(x,k?) presents the oscilla- 
tions near x c when /c 2 — ► Notice that these oscillations are random since 
the evolution is not ordered on the transverse momentum space. For illustrat- 
ing the details of the oscillations, we re-plot a curve with x = 7.6 x 10~ 7 of Fig. 
16 in Fig. 17a, where linear scale replaces logarithmic scale on the longitudinal 
axis. 

(hi) The sudden change of a solution is an interesting phenomenon in non- 
linear evolution system, in particulary, this behavior perhaps relates to chaos. 
We study the dependence of the solution on the input conditions. For this 
sake, we compute a similar solution as in Fig. 17a but the starting point is 
moved from xq = 0.001 to x = 0.00105. The results in Fig. 17b show that 
the oscillation structure of F(x,k 2 )/h^ ~ A; 2 is sensitive to the starting point 
of the evolution, although the global behaviors of the curves are similar. 

(iv) We change the parameter 40 to 45 and 50 in Eq. (4.7) and compare 
these results in Figs. 18a, 18b and 18c, respectively. One can find the obvious 
difference in the oscillation structure if using different input parameters. 

(v) The improvement of the precision in the computation may aggravate the 
chaotic oscillations since the increasing samples perturb the distributions at 
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every step in the evolution. In contrast to it, if the above mentioned oscillations 
are arise from the calculation errors, such oscillations will be suppressed with 
the increasing precision. In Fig. 19 we compare a curve with x = 7.6 x 1CT 7 
in Fig. 16 (see Fig. 19a) with a same solution but using double calculating 
precision (see Fig. 19b). One can find that the oscillations are aggravated 
with increasing precision. On the other hand, some global properties of the 
solution, for example, the position of x c is unchanged. 

(vi) The above mentioned properties of the MD-BFKL equation show a 
typical chaotic behavior: the distribution F(x, fc 2 ) ~ k 2 presents the random 
oscillation in a certain kinematic domain, and its structure sensitively depends 
on the input conditions. Thus, we suggest that the MD-BFKL equation has a 
chaotic solution at small x. 

(vii) Now let us point out that the above mentioned chaotic solution in 
the MD-BFKL equation directly arises the strong shadowing. We calculate 
the contributions of a nonlinear regularized factor in Eq. (4.5) (i.e., the factor 
inside the large brackets) to AF(x : k 2 ). The results are shown in Fig. 20, 
where we define B as the contributions from the part of 

81 aj r°° dg j k 2 F 2 (x,k! 2 ) _ k! 2 F 2 (x,k 2 ) } 
4 nR% 4 2 , it' 2 I it' 2 W 2 - k 2 1 £' 2 1 U 2 - k 2 1 J { } 

to AF(x, k 2 ) via Eq. (4.8). The first and second terms in Eq. (4.8) are 
the contributions from the real and virtual processes, respectively. Figure 20 
shows that the contributions of two terms are almost cancelled each other and 
remains the smaller positive contributions in the region of x ^> x c = 7.2 x 10~ 7 , 
where F(x,k?) is a smooth decreasing function of ic 2 when fc 2 — > k%. However, 
the chaotic solution of Eq. (4.5) near x ~ x c changes the situation: the factor 
{...} in Eq. (4.8) becomes ~ F(x : k/ 2 )/k/ 2 . This leads to a larger negative 
contribution of Eq. (4.8) to AF(x,fc 2 ) in the region of the violent oscillation 
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at x ~ x c and A;' 2 — >■ A; 2 ,, i.e., form the strong shadowing. It seems that a similar 
enhancement also should appear in the positive linear regularized factor in Eq. 
(4.5). However, because the oscillation of F 2 (x, fc 2 ) is stronger than that of 
F(x,k 2 ), the evolution stopped before the linear regularized factor enhances. 
Therefore, we conclude that the gluon disappearances in the solution of the 
MD-BFKL equation is caused by a strong shadowing effect, which exists in 
the chaotic solution of this evolution equation. 

The above mentioned calculations of the MD-BFKL equation contain some 
undetermined parameters. They are the correlation parameter in Eq. (4.5) 
and the normalization constant (5 in Eq. (4.7). However, we find that these 
uncertainties don't change either the chaotic properties of the solutions or the 
conclusion about the gluon disappearance, although they hinder us to predict 
the value of x c . 

Other relating problem is the value of k^. In the case of the running a s , the 
value of kl = O.OlGeV 2 goes into the nonperturbative region of k 2 , where the 
application of the perturbative QCD is questionable. Usually, we can introduce 
a lager cut-off A^ to avoid the infrared region. Our calculations show that the 
chaotic solution of the MD-BFKL equation always occurs near the value of 
k% even it becomes large. For example, if k% = lGeV 2 and (3 = 1, a similar 
oscillation in F(x,k 2 ) will appear near x c ~ 3.8 x 1CT 7 and k_l ~ lGeV 2 . 

It is interesting that according to the dipole picture of the BK equation, 
the nonlinear term in this equation can be explained as the double scattering 
of the probe on two connecting sub-dipoles (see Fig. 21a). One can image 
that the probes also can double scatter on a same sub-dipole, for example, 
as shown in Figs. 21b and 21c. These processes shall add a term to the BK 
equation, which is proportional to ~ N 2 (x,x). In this case, the contributions 
of the corresponding virtual processes are necessary and the same effect in the 
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MD-BFKL equation also possibly exists in this kind of the BK equation. 

Let us briefly discuss the measurements of the new effect of the MD-BFKL 
equation. By the look of the probe in a DIS process, all partons have the 
transverse size ~ l/Q, where Q is the probe virtuality. According to the 
pictures in Fig. 1, we schematically draw the kinematic regions of four evo- 
lution equations in Fig. 22, where the transverse direction points the linear 
approximation, while the longitudinal axis spans the kx- and collinear factor- 
ization schemes. The BFKL region is below than the DGLAP region since the 
cold spots are frequently formed in the lower Q- and smaller x-region. Note 
that the evolution dynamics in the overlapping region of the MD-DGLAP and 
BFKL equations is the BK-like equation, which can be structured by the ker- 
nels of these two equations [17]. One can find that the MD-BFKL equation 
works following the BK-like equation, where the gluon distribution has reached 
saturation. Thus, we expect that the gluon disappearance occurs below the 
saturation behavior. 

We focus on the region Q 2 < Q 2 s {x) (Q 2 s (x) ~ a few GeV 2 at x < 10~ 3 
and Q s {x) is called the saturation scale [7]), where the MD-BFKL equation 
dominates the evolution procuresses. We assume that the starting point x in 
Eq. (4.7) is irrelevant to the values of k 2 . Using Eq. (2.41) and the results in 
Fig. 13, we conclude that the (integrated) gluon distribution has the form 

G(x, Q 2 < Q 2 S ) = at x< x c . (4.9) 

Predictions for the radiative quark distributions are then made by convo- 
luting the unintegrated gluon distribution with the photon-gluon fusion quark 
box using 

x(q(x,Q 2 <Q 2 ) + q(x,Q 2 <QD) 
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f 1 , fQ 2 dk 2 [*? dk' 2 „. x ,«,a Sn . N 
= 2 .W "f/ Zi 2^ ' (410) 

which predicts that the sea quark distributions gradually approach to zero with 
x — > 0. It is interest that the valence-like gluon and sea quark distributions 
have been used in some input parton distributions in the QCD-global analysis 
of the present experimental data [19]. 

On the other hand, at the region Q 2 > Q 2 (x), the DGLAP or MD-DGLAP 
equations will replace the MD-BFKL equation. One can find that the sunken 
part in G(x < x c , Q 2 ) at Q 2 3> Q 2 will be fast filled up by the soft gluons via 
the DGLAP dynamics (A-9). 

In our above numerical analysis, the 4-gluon correlation function (x, k) 
is assumed as the square of the unintegrated gluon distribution since the con- 
dition Eq. (3.31). Of cause this restriction condition is approximal due to Eq. 
(3.21). A smaller deviation from Eq. (3.31) is possible. For this sake, we take 
the recombination probability as an exponential derivation form, i.e., 

R F {k 2 ^ d )=^Be-p k2 - k ^\ (4.11) 
where B is the normalization constant 

B=°[2-e-°]-\ (4.12) 

and a is a lager parameter implying the smaller deviation from Eq. (3.31). 
Note that R F (tf,l£d) oc 5(^ d - k 2 ) if a > 1. We recalculate the MD-BFKL 
equation using Eqs. (4.7) and (4.11), i.e., 

dF(x,k 2 ) 
dx 

3a s k 2 r°° dti 2 J F{x,ti 2 ) - F(x,k 2 ) F(x,k 2 ) 1 

~~M 1^1 W 2 -k 2 \ + ^ + ip J 
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l F(x,k 2 )F(x,£ d ) 



41 dkj d B 2 e 



-fj\k 2 -k 2 ad \ 



k' 2 \k' 2 -k 2 \ 

2\ cV™ k2 



F(x,k 2 )F(x,ki d ) 



^k 4 + 4k! 4 



(4.13) 



where 



51 = F [2-eT, (4-14) 

5 2 = p[2-e- a ]- 1 - (4-15) 

For example, we take a — 1, 5, 20 and oo, respectively. The results for k 2 = 
5 GeV 2 are shown by the solid curves in Fig. 23. The dashed and point curves 
are the BFKL- and BK-like-solutions, respectively. We find that the drooping 
position x c of the gluon distribution moves toward the smaller value of x when 
the parameter a in Eq. (4.11) decreases, i.e., the deviation from Eq. (3.31) 
becomes larger. We understand that the change of the nonlinear terms due to 
^ ad 7^ k 2 reduces the chaotic solution of the MD-BFKL equation. Although 
we have emphasized that a lager deviation from Eq. (3.31) is unreasonable, 
the above mentioned results increase the uncertainty in the predictions of the 
position, where the gluons disappear. We also draw the solutions of the BK 
equation (point curve) with the same input distribution. One can find that 
before the gluon distribution dropping, the screening effects in the MD-BFKL 
equation is stronger then that in the BK-like equation. 

Obviously, the MD-BFKL equation is derived in a partonic framework, we 
regard the present form of this equation is a possible evolution model rather 
then a proper QCD dynamics. An available QCD evolution equation should 
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be structured based on the renormalization group equations, where the Wil- 
son coefficients are separated from the anomalous dimensions. Although the 
equivalence between the TOPT and renormalization group methods in the 
leading twist DGLAP equation has been justified, the evolution kernels of the 
MD-BFKL equation should be proofed. On the other hand, we have not relate 
to the factorization schemes including the higher twist operators, this difficult 
problem has gone beyond a naive evolution model. As we well know that the 
next-to-leading logarithmic (NLL) corrections to the BFKL equation lead to 
the large numerical effects. A particularly worrying consequence is that the 
eigenvalue becomes negative and this would indicate an unphysical negative 
cross section [20]. Many studies have been devoted to the analysis of the BFKL 
approach at NLL accuracy in order to obtain reliable predictions [21]. An im- 
portant question is whether the MD-BFKL equation is repeating a similar 
mistake as in the above mentioned NLL BFKL equation? We think there is 
an essential difference between two evolution equations: the large corrections 
in the linear NLL BFKL equation can be improved using the reorganization of 
the perturbation expansion, while the chaotic solution in the MD-BFKL equa- 
tion is a peculiar phenomenon of the nonlinear dynamics and it is observable. 
Of cause, a complete MD-BFKL equation should include the NLL corrections 
with running coupling. This is our next aim. 

In summary, an evolution equation incorporating the corrections of the 
gluon fusion to the BFKL equation is proposed in a unified partonic framework. 
This modified BFKL equation predicts a stronger shadowing, which rapidly 
suppresses the gluon density and even leads to the gluon disappearance below 
the saturation region. We find that this unexpected effect is caused by a 
possible chaotic solution in the new equation. 

Acknowledgments: We thank Z.H. Liu and J. Zhou for useful discussions in 
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Appendix A: The TOPT-cutting rules 

The TOPT-cutting rules state the possible relations among the different 
cut diagrams in DIS, which are developed in our previous works [9,12,16]. The 
idea of the TOPT-cutting rules are follows. The propagators in a DIS process 
contribute either their forward or backward components in the Bjorken limit, 
i.e., both them are the mass-shell. Therefore, one can replace a propagator by 
a pair of in- and out-lines, thus, the propagator can be cut. 

For illustrating the TOPT-cutting rules, we write a complete derivation of 
the DGLAP equation for gluon. The contributions of the real diagrams (Fig. 
24a) to the hard part of the evolution equation are 

H(probe*l — > -probe* I) 

1 „ 1 1 1 1 rf3/ ' 
-Mi^uv M,* hV 

O 771 771 TP TP o TP TP TP TP (— >/Ci 



2E l 2E k E L -E k - E v 2E k E x -E k - E v (2ir) 3 2E v 

1 d 3 k' 
x— \M probe * k ^ k ,\ 2 (27r) 4 S(p probe * + k-k') 



• AE probe ^ prooe 1 ' ^' ooe ' (2vf2E kl 

dl 2 

= IC DG lap(xi, x 2 , x 3 , a s )dxs—^5(x 2 - x B )dx 2 C (probe* k — > k'), (A - 1) 

where the factor 1/(2E{) is extracted from the definition of the gluon distribu- 
tion g(Q 2 , xb), and 



JC DGLA p(x 1 ,x 2 ,x 3 ,a s )dx 3 



Eh , ,r> 

E 



Ei — E k — Eii 



dl 2 

T 2 

l i 2 dH> 



2E h 



(27r) 3 2iV 



(A -2) 



C (probe* k — > k')5(x 2 — xs)dx 2 



M probe * k ^ k t\ (2n) 5(p probe * + k - k') , (A - 3) 



o£j k hj probe * ynx)°ihj k i 
is the contributions from da(probe*k — > A;'). Thus, 
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f dl 2 

Ag(x B ,Q 2 ) = J -p L dx 1 dx 2 dx 3 JC DGLAP (x 1 , x 2 ,x 3 , a s )g(x 1 ,l 2 L )5(x 2 -x B )5(x 1 -x 2 - 

(A-4) 

where we inset 8{x\ — x 2 — x 3 )dx\. We obtain 

dg(x B ,Q 2 ) 
d\nQ 2 

= j dx 1 dx 2 dx 3 dzfC DGLA p(x u x 2 ,x 3 ,a s )g(x 1 ,Q 2 )5(x 2 - x B )S(x 1 - x 2 - x 3 ) 
= [ —ICDGLApf—^s^gix^Q 2 ). (-4-5) 

J xb X\ \ X\ J 

On the other hand, the contributions from one of the virtual diagrams, for 
say, Fig. (23b) are 

H(probe*l — > probe*l) 

11,, 11 1 „ 1 1 d 3 l' 

-M,^ kl , — — - —M kV . 



2 2E t 2E k 2E V E L — E k — E v 2E t E k + E v — E\ (2ir] 

1 d 3 k' 
x — \M probe * l _ tk ,\ 2 (2ir) 4 5(p probe * + l-k') 



3 



1 ^ , s . dl 2 

where 



= --K. DGLA p(x 1 ,x 2 ,x 3 ,a s )dx 3 - r 5(x 1 - x B )dx 1 C(probe*l -> ft;'), (A - 6) 



C (probe* I — > — x B )dx\ 

1 d 3 ft;' 

'M pr . of)e */^| 2 (27r) 4 (5(pp ro f )e . + / - ft/) — , (A - 7) 



oHillliprobe* \4Tf) *£>k' 

is the contributions from dcr(probe*l — > ft;'). Therefore, we have 

<91nQ 2 

= ~\J dx 1 dx 2 dx 3 K, DGLAP (x 1 ,x 2 ,x 3 ,a s )g(x 1 ,Q 2 )5(x 1 -x B )5(x 1 - x 2 - x 3 ) 
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1 f 1 

-7,9(xb,Q 2 ) / dzK, DGLAP (z,a s ). 

2 Jo 



(A -8) 



The complete evolution equation for the gluons is 



dg(x B ,Q 2 
dlnQ 2 



where 



/ —fc DGLA p(—,a}\g{xi,Q 2 )-\g(x B ,Q 2 )( dzJC DGLAP (z, a s ), 

Jx B X 1 \Xl J 2 Jo 

(A -9) 



DGLAP 



a,N n 



71 



1-Z 



+ 



1 



+ z(l-z) 



(A - 10) 



Note that the factor 1/2 in Eq. (A-8) is the considerations of the symmetry 
under exchange of two internal gluons in the virtual diagrams (see Ref. [9,12]). 
The expressions of the evolution equation in the TOPT form, without the 
calculations of the matrixes, show that the real and virtual diagrams contribute 
the same evolution kernel but with the different factors. This simple form of the 
equation consists with a more complicated derivation of the DGLAP equation 
using the covariant perturbation theory in [11]. In fact, from Eqs. (A-9) and 
(A-10) we have 



dx B g(x B ,Q 2 ) a s N c 



dlnQ 2 



71 Jx B 



dz 



z 1 - z . 
+ +z(l-z) 



1-z 



z 



Xig(xi,Q 2 ) 



la ' Nc -x B g{x B ,Q 2 ) fdz 



2 7T 







z 1 — z . 
+ + z(l-z) 



1-z 



(A -11) 



where z in the first term of the r.h.s is read qsxb/xi- Note that Jq 1 dzz/(l—z) 
Jq 1 dz(l — z)/ ' z, we obtain 



dx B g(x B ,Q 2 ) _ a s N c f 1 
dlnQ 2 

a s N c 



71 Jx B 



dz 



zx 1 g{x 1 ,Q 2 )-x B g{x B ,Q 2 ) {1 - z){l + z 2 ) 2 



1-z 



—x B g(x B ,Q 2 ) [ dz 

71 JO 



V -z(l - z) 

1-z 2 v ; 



a a N, 



z 

71 Jx B 1 — 
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a a Nc f 1 

71 Jx B 



dz 



X!g(xi,Q 2 )z - x B g(x B ,Q 2 ) (1 - z)(l + z 2 ) 2 
1 x^x^Q ) 



g s N c 

7T 



1-Z 
11 



12 



+ ln(l - x) 



x B g(x B ,Q 



(A - 12) 



This equation is equivalent to the following traditional form of the DGLAP 
equation 



dg(x B , Q 2 
d\nQ 2 



1 ^/c^r (—,«.) g{x x ,<?)), 

x B X\ \X\ 



(A - 13) 



where 



^traditional ( 



'DGLAP 



(*) = 



a«N r . 



7Y 



_) +__ + «!-,) + ^(1-,) 



(.4 - 14) 



In fact, using 



o (1 - z)„ 



1-2 



(A - 15) 



Equation (A-14) becomes 



dx B g(x B ,Q 2 ) a s N c r 1 



dlnQ 2 



7T ii B 



1-Z 



1-2 



+ 2(1-2) 



11 

H a 

127T ' 



<N c x B g(x B ,Q 2 ) - ^-^5(1 - 2) / d2-^ 

7T JO 1 — 



« S A C /"I 



7T 



G?2 



x 1 g(x 1 ,Q 2 )z - x B g( x B ,Q 2 ) | (l-2)(l + 2 2 ) 2 

t Xigyxi, j 



+ 



a,A n 



7T 



1-2 

— + \n(l-x B ) 



x B g(x B ,Q ) 



(A - 16) 



One can find that IR divergences of the real diagram at the end point (x B — > 1) 
are exactly cancelled by the the contributions of the virtual diagram. 

We emphasize that an important property is used in the above mentioned 
derivation: two propagators in Fig. 24b only contribute their TOPT-forward 
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components in the Bjorken frame [12]. Thus, the propagator on mass-shell 
can be replaced by a pair of in- and out-lines, i.e, the propagator can be cut. 
A similar way was generalized to include the multi-initial partonic processes, 
in which the resummations of more complicated cut diagrams are included 
[9,12]. Besides, in the nonperturbative correlation function, the propagator on 
the light cone coordinator y~, which connects with the perturbative scattering 
amplitude either only keeps either its backward component (in the collinear 
factorization scheme [16]) or forward component (in the W — W approxima- 
tion), i.e., this propagator is on mass-shell and can be cut. Summary, the 
TOPT-cutting rules can be summarized as follows: 

Assuming a general structure of the cut diagrams G(N) in TOPT is 



g,w= n v — /-v F n v f X — w- 

left-vertices ^a'ei+1 ^a' l^aEi right-vertices <^bej l^b'ej+1 

x n *(e*/=°) n 2^ 

vertices-k \ f ) states— c c 

where N G is an overall numerator-and-symmetry factor and is independent of 
the cut 7; i and i + 1 (or j and j + 1) are the time-ordered lines on the left- (or 
right-) vertices; S(J2f x f — 0) is the conservation of longitudinal momentum 
at the vertex. Thus, the contributions of the cut diagrams in the sum X) 7 C 7 
have the identical integral kernel with only the following different factors R: 

R = (±) x (1, 1) x 5(xp - x B ), {A - 18) 

(a)The sign in the first factor of (A-18) is determined by the energy deficits 
in (A- 17). 
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(b) The second factor takes a value of 1/2 if the probe- vertex inserts in the 
initial line. 

(c) 5-function relates to the probe vertex. 

(d) When the cut line moves its position, the contributions of the final 
states in (A-17) will change the momentum-symbols, but dont change the 
structure of the intermediate state. 

(e) The cut line can also cut the nonperturbative matrix elements with 
multi- initial partons: the initial parton line on the light-cone can be moved 
from one-side of the cut to another side. 



46 



Appendix B: A unity of the QCD evolution equations 

Through the derivations in this work, a set of the QCD evolution equations 
at small x are presented in a unified partonic framework. We summary them 
as follows. 

The DGLAP equation (2.41) 



2 dG(x,Q 2 ) a s N c f 1 dx 1 

71 Jx Xi 



dQ 2 



[ ^lG(x u Q 2 ), (B-la) 

Jx Xi 



or 



AG(x,Q 2 ) 
(see Fig. 25a); 



7T Jx X 1 



a s N c r 1 dx x fQ 2 dk 2 



k 



-G{ Xl ,k 2 ) 



The BFKL equation (2.39) 



(B - lb) 



-x- 



dF(x,k 



aO) 



dx 



a s N r rd 2 ^ 



a s i\ c r u k 

IT J IT 

or its real part 



k 2 k 2 

^ab^Gb 



(B - 2a) 



a s N c f 1 dx x f d 2 ka b , 

2 , 2 r \ X l->!±ab) 



7T Jx X\ J TT 

(see Fig. 25b); 

The MD-DGLAP equation (3.36) 



k A k z 



(B - 2b) 



dGjx, q< 
dlnQ 2 
a s N c f 1 dxi 

Xi 



S N C f 1 dx x 2 

7T Jx Xi 



18a 2 N 2 /-V2 dx x 



ttRIQ 2 N 2 



c —\ —G 2 (x l7 Q 2 ) 

— 1 Jx/2 X\ 
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or 



3Qai 



ttR 2 n Q 2 N 2 



N? f 

j 2 - 1 Jx 



V2 dXi 
Xi 



G 2 (x 1? g 2 ), 



(B - 3a) 



AG{x,Q 2 ) 

<\„X, f 1 dxi fQ 2 dk 2 
~k 



IT Jx X-i J 



2 G(x 1 ,k 2 ) 



18a 2 N 2 /-V2 dx x fQ 2 dk 2 1 



ttR 2 n Q 2 N 2 - 1 Jx/2 Xl J k 2 k 2 
36a 2 N 2 /-V2 d Xl fQ 2 dk 2 1 



r 1 / z dxi r 

Jx Xi J 



itr 2 n q 2 n 2 -iJx xi j e k 2 

(see Fig. 25c); 

The MD-BFKL equation (3.46) 



G (xi,k ) 



(B - 36) 



= a s N c f d 2&h 

IT J IT 



-X- 

kbc kbo 



dF(x,k m ) 



k 2 k 2 



dx 

F(x,kbc) - ^Ffakw) 



ISa 2 N 2 
ttR 2 n N 2 - 1 



d 2 kjbc kbo 
7T ki,„k 



36a 2 N 2 
~*R%N 2 - 

or its real part 



— bc^-cO L— be 



-irl" 2 i-.k ) - 



kV V2'- fec ; 2k 2 



d 2 k k 2 
71- p P 



z I. 

2 



T2-F 2 (2U 6c ) 



2/V 2 



(B - 4a) 



AF(x,/c w ) 



_ a s N c r 1 dx± r £% c /c 2 

IT Jx Xi J IT 



18a 2 Nl 
ttR 2 n N 2 - 



36a 2 N 2 
ttR 2 n N 2 - 1 

(see Fig. 25d). 



V2 aXi 
/2 X! 

V 2 dxi f d 2 kb c fcj 



,2 t 2 F{ X likbc) 



r L i z dx 1 r a jz 

1 Jx/2 Xi J TT 



rl 2 h h 2 



1 



/V 2 /V 2 /V 2 
HcHa He 



F 2 (x 1 ,k bc ) 



IT kbckcQ hbc 



{B - 46) 
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One can find following interesting relations among these equations: 
(1) The DGLAP and BFKL equations have the same evolution dynamics 
(i.e., the gluon splitting), which evolve on the Q 2 -axis and the (& 60 , /c b0 )-plane, 
respectively. This leads to following relationships between the real parts of 
Eqs. (B-lb) and (B-2b): 



dk (Pkgb k a0 
k 77 k nh kn 



(5-5), 



G(x,k 2 ) F(x,k ab ). (B-6) 

The MD-DGLAP and MD-BFKL equations also have the same evolution 
dynamics (i.e., the gluon recombination), which evolve on the Q 2 -axis and the 
{kbo, iy>o)~pl ane > respectively. Therefor, one can expect that there are similar 
relationships like Eqs. (B-5) and (B-6). In fact, comparing Eqs. (B-3b) with 
(B-4b), we have 



dfc d —bo 



k 7T k bc k c0 



(B-7) 



G(x,k 2 ) F(x,k bc ), (B-8) 
and an extra relation for the power suppression factor 



1 

72 N ' 7 2 7.2 

K Kbc IHo 



(B-9) 



Besides, IR divergences both in the evolution kernels of the BFKL and MD- 
BFKL equations are cancelled by using the summations of the contributions 
from the real and virtual processes. We emphasize that the same summation 
is also used in the DGLAP equation (see Eq. (A-9)) and the MD-DGLAP 
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equation [9,12] but the virtual contributions are neglected at small x in Eq. 
(B-la) and they are cancelled due to the symmetry in Eq. (B-3a), respectively. 

(2) From Fig. 1 one can find that if the interferant amplitudes are negligible 
the BFKL and MD-BFKL equations reduce to the GDLAP and MD-DGLAP 
equations, respectively. 

(3) Both the MD-DGLAP and MD-BFKL equations contain the correc- 
tions from the interferant and un-interferant processes and they contribute 
the shadowing and antishadowing effects, respectively. 

(4) All equations have the smooth solutions except that the MD-BFKL 
equation contains the chaotic solution, which relates to the following special 
structures: the MD-BFKL equation contains the nonlinear distributions com- 
bining the regularized BFKL-kernel. 

The above mentioned discussions show that the MD-BFKL equation (3.36) 
consists with the DGLAP, BFKL and MD-DGLAP equations. In other words, 
the present form of the MD-BFKL equation is reasonably determined by the 
DGLAP, BFKL and MD-DGLAP equations. 

On the other hand, the BK equation for the unintegrated gluon distribution 
reads [22] 



Obviously, it is difficult to combine this form of the BK equation with the 
above mentioned united evolution equations. 




dF(x,k 2 ) 
dx 

3a s k 2 r°° dg 2 \ F(x, ti 2 ) - F(x, k 2 ) F(x, k 2 ) 




(B - 10) 
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Figure Captions 

Fig. 1 The corrections of the initial gluons to the basic amplitude of the 
DGLAP equation and they lead to (b) BFKL-, (c) modified DGLAP- and (d) 
modified BFKL equations, respectively. 

Fig. 2 The TOPT-diagrams consisted by the elemental amplitudes in Fig. 
lb. These diagrams lead to the real part of the BFKL equation. The dashed 
lines are the time ordered lines in the TOPT. 

Fig. 3 The QCD correlations between two initial gluons. 

Fig. 4 A sub-vertex in a complex Feynman diagram. 

Fig. 5 The Virtual diagrams corresponding to Fig. 2. 

Fig. 6 A cutting diagram originating from Fig. Id. 

Fig. 7 The TOPT-diagrams consisted by the elemental amplitudes in Fig. 

Id. 

Fig. 8 The Virtual diagrams corresponding to Fig. 7. 

Fig. 9 One of the TOPT-diagrams for the interference processes, which 
have the same order as Fig. 7. 

Fig. 10 The Virtual diagrams corresponding to Fig. 9, they contain a 
similar evolution kernel but with a different numerical factor according the 
TOPT-cutting rule. 

Fig. 11 The TOPT-cutting rules shows a simple relations among virtual 
diagrams in Figs. 7 and 9. Thus, all diagrams in Figs. 7, 8, 9 and 10 contain 
a similar evolution kernel but differ from the different numerical factors. 

Fig. 12 The model for the multi-gluons correlating function based on the 
TOPT-cutting rules. The propagator inside the cold spot is forward and on 
mass-shell at the W — W approximation, while the correlations to the cold 
spot from the other part of the unperturbative matrix are neglected. Thus, 
the correlating function can be cut. 
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Fig. 13 x-dependence of the unintegrated gluon distribution in the MD- 
BFKL equation (4.5) for two different values of k 2 . The dashed curves are the 
corresponding solutions of the BFKL equation. The results show that F(x : k?) 
suddenly drops near a critical value of x ~ x c ~ 7.2 x 1CT 7 . 

Fig. 14 A same solution with k 2 = 5GeV 2 in Fig. 13, but where the 
contributions from the nonlinear shadowing part in Eq. (4.5) are separately 
indicated. The results show that the shadowing effect increases suddenly near 
x c . 

Fig. 15 Similar to Fig. 14 but for the solution of Eq. (4.6). The results 
show that the shadowing effect grows gradually. 

Fig. 16 /^-dependence of the unintegrated gluon distribution (scaled by k 2 ) 
for different values of x. The dashed curves are the corresponding solutions of 
the BFKL equation. The results present the oscillations near x c when A; 2 — > k^. 

Fig. 17 Comparisons of the /^-dependence of the unintegrated gluon dis- 
tributions at x = 7.6 x 1CT 7 in the MD-BFKL equation (4.5) but evolving from 
two different stating point x . The results show that the oscillation structure 
is sensitive to the starting point x of the evolution. 

Fig. 18 Comparisons of the /^-dependence of the unintegrated gluon distri- 
bution at x = 7.6 x 1CT 7 in the MD-BFKL equation (4.5) but using the input 
(4.7), where the parameter 40 (Fig. (a)) is replaced by 45 (Fig. (b)) and 50 
(Fig. (c)), respectively. One can find the obvious difference in the oscillation 
structure if using different input parameters. 

Fig. 19 (a) The curve with x = 7.6 x 10~ 7 in Fig. 16; (b) the same as (a) but 
is computed by the double precisions. The results show that the oscillations 
are aggravated with increasing precision. 

Fig. 20 The contributions of the nonlinear regularized factor in the MD- 
BFKL equation (defined in Eq. (4.8)) to AF{x,k 2 ). 
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Fig. 21 (a) The dipole model of the BK equation, where probes are scat- 
tered off the two dipoles; (b) and (c) The probes are scattered off the same 
dipole, which are neglected in the BK equation. 

Fig. 22 Schematic kinematic regions of four evolution equations, where the 
transverse direction points the linear approximation, while the longitudinal 
axis spans the kx- and collinear factorization schemes. The evolution dynamics 
in the overlapping region of the MD-DGLAP and BFKL equations is the BK- 
like equation. 

Fig. 23 The solutions of Eq. (4.13) for A; 2 = 5 GeV 2 , where the correlations 
of two different initial cold spots in the MD-BFKL equation are considered: the 
solid curves from left to right are a — 1, 5, 20 and oo. The dashed and point 
curves are the corresponding solutions of the BFKL and BK-like equations, 
respectively. 

Fig. 24 TOPT diagrams for the DGLAP equation: (a) real diagram, (b) 
and (c) virtual diagrams. 

Fig. 25 A part of the amplitudes for four evolution equations based on the 
model of Fig. 12. 
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Fig. 9 
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